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1. Introduction 



The six vertex model, a well-known integrable model of two-dimensional statistical 
mechanics, has been solved with various types of boundary conditions: in the original solu- 
tions, with periodic boundary conditions (PBC) then more recently with anti-periodic 
boundary conditions and "domain wall" boundary conditions (DWBC) ||,|5],|6|,[7],|8[] . In- 
terestingly enough, in the latter case, even "bulk" quantities turned out to be different 
than with PBC. One possible way to understand this is to notice that the six-vertex model 
possesses at each vertex a constraint: equality of the number of incoming and outgoing ar- 
rows. This conservation of arrows renders the usual theorems of statistical mechanics which 
ensure the existence of a thermodynamic limit with bulk quantities that are independent 
of the boundary conditions inapplicable (some energies being infinite); and it effectively 
creates a non-locality of the degrees of freedom (or of the moves from an algorithmic point 
of view) which is expected to create sensitivity on boundary conditions. This raises the 
general problem of the effect of boundary conditions on the thermodynamic limit of the 
six-vertex model. Important questions such as the computation of the bulk free energy, as 
well as of some local quantities (local polarization) will be addressed here. We shall try 
to give a qualitative understanding of the physical phenomena involved, as well as some 
explicit calculations whenever they are possible. 

We shall consider in this paper fixed boundary conditions (FBC) from which most 
boundary conditions can be derived. As we shall see, the main idea is that these boundary 
conditions induce local polarizations.0 These are responsible for spatial phase separation: 
considering a continuum limit with a proper scaling for the domain in which the model is 
defined and for the limiting conditions at its boundary, we shall see that several phases 
may coexist in this domain. Depending on the phase diagram (that is on the value of the 
Boltzman weights defining the model), the phases have definite boundaries in the case of 
second order phase transition, whereas in the case of first order the phases freely mix. 

The techniques we shall use are strongly inspired by the recent developments in the 
field of dimers / domino tilings |9|,10,ll JT^1 which are themselves related to general ideas 



in random tilings |13[]. As was pointed out in this context in J7]||], domino tilings are 
a particular case of the six-vertex model, corresponding to a special set of values of the 
Boltzman weights. 



1 FBC fix in particular the total polarization, but the statement is stronger. 
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The plan of the article is as follows. In section 2, we shall define the six-vertex model 
with general FBC and try to justify a conjectured variational principle for the description 
of its thermodynamic limit. Then we shall analyze the phase diagram of the model in the 
three regimes that it possesses and give explicit equations for the variational problem in 
some cases (section 3) . Finally, we shall conclude in section 4 with some comments on the 
applications of this variational principle (in particular in relation to conjectures made in 
0H)) an d some open questions. 



2. The six-vertex model and its variational principle 

The six- vertex model is defined on a regular square lattice. For fixed boundary condi- 
tions (FBC) the lattice will be contained in a domain T> of the plane which is assumed, for 
simplicity, to be convex. For technical reasons which will become clear below, we must also 
consider in parallel the model with periodic boundary conditions (PBC), in which case the 
lattice is on a torus. The lattice spacing is called S. 

The configurations of the model are obtained by assigning arrows to each edge of the 
lattice (see Fig. 1). The FBC mean that the arrows at the boundary of T> are supposed to 
be fixed once and for all. 





>- 


/ \ 


,\ 
,<r~ 
/ \ 


7 / 
< 

/ / 


< 






■A 


\ \ 


> 
/ \ 


> 
1 \ 


> 

/ \ 


< 

/ \ 


/ / 




>- 


-<r 

/ \ 


1 \ 


/ \ 


J \ 


/ / 


\ / 


\ : 


< 


< 
/ \ 


< 
1 / 


> 
\ \ 


> 
/ \ 


> 
/ / 


> 
\ / 


V 


< 

\ 


< 
/ \ 


< 
1 / 


< 
\ / 


> 

\ \ 
< 


> 
/ \ 


/ 




>- 
"a 


> 
/ \ 


> 
/ \ 


< 
/ / 


\ / 


> 
\ \ 





Fig. 1: A configuration of the six- vertex model. 
The partition function is then obtained by summing over all possible configurations: 

z= yi n w ^ (2- 1 ) 

arrow configurations vertex veD 
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where the statistical weights w(v) are assigned to each vertex v of the lattice. These are 
given in terms of the arrows around the vertex according to 
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(2.2) 
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where a, 6, c are positive real numbers. These being the only configurations that respect 
the conservation of arrows, all the other weights are zero. 

One can consider the arrows as vectors in the plane with unit length; then the polar- 
ization of a subdomain V C V is defined as 



v ' edge eel?' 



where Af(V) is the number of vertices in T>', and p(e) is the arrow of edge e. Note that 
with FBC, the total polarization per vertex P = P(T>) is fixed, but not with PBC. In the 
latter case it is natural to introduce an electric field E coupled to the polarization: 



z(e)= ex p(^- e jfc)) n 

arrow configurations ed s e eeV 



w(v) 



(2.3) 



vertex vET> 



(we set for now the temperature to one. it will be restored when needed). In the case of 
PBC it is known that the thermodynamic limit is well-defined; that is, independently of 
the shape of the torus, when one sends its linear sizes to infinity there is a unique free 
energy per vertex 

F(E) = -km o ±\ogZ(E) (2.4) 

where Af = Af(V) is the number of vertices of the lattice. One also defines for future use 
the Legendre transform of F: 



G(P) = F(E) + E-P 
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(2.5) 



— * — * 

where P and E are related 




(2.6) 



As is well-known, G is the free energy per vertex at fixed total polarization. 
In the case of FBC, one also defines the free energy per vertex to be 



Gfbc = — 77 log Zybc 



(2.7) 



where the subscript FBC reminds us that G depends on a choice of the arrows at the 
boundary of V. 

Finally, we need to define the height function h associated to a configuration of the 
six-vertex model with FBC; it is a function defined on the faces of the lattice, such that 
when one moves from a face to one of its neighbors, h increases by 5 if the arrow in between 
points right (and of course decreases by 5 if it points left), h is only defined up to a global 
constant. Similarly, there is a boundary height function ho defined as above on the faces 
surrounding T>, and which is fixed up to a constant when the boundary arrows are fixed. 
In the end, the whole model can be recast as a height model: the configurations are now 
aL- valued functions h defined on the faces of the lattice, which satisfy \hf — hf\ = 8 for 
neighboring faces / and /', and the boundary condition h\gx> = ho- The Boltzman weights 
are the direct translation of those of the original model. 

Now we are ready to take the continuum limit for FBC. It is obtained by sending 
the lattice spacing S to zero, keeping the domain V fixed. What happens to the boundary 
conditions? Since they can be specified by the boundary height function ho, we assume that 
the latter converges to a given function ho defined on the boundary dV. It is important to 
notice that typical configurations of the height function will also be described by smooth 
functions defined on V. In fact the constraint \hf — hp\ =6 guarantees that any limiting 
height function will satisfy the 1-Lipschitz condition \h(x, y) — h(x', y')\ < \x — x'\ + \y — y'\, 
where x and y are coordinates along the directions of the lattice, so that h will be (almost 
everywhere) differentiate. 

More precisely, let us consider the average height function (h(x,y)), which plays the 
role of one-point correlation function in the model. In the continuum limit, (h) becomes a 
macroscopic quantity. In fact, its derivative is nothing but the average local polarization: 



d_ 

dx 
d_ 

dy 



(h(x,y)) = {P y (x,y)) 



(h(x,y)) 



< 



P x (x,y)) 



(2.8a) 



(2.86) 
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P(x,y) being denned by 




,y) 



where T>' is a macroscopic patch around (x,y) whose size tends to zero. 

The important consequence of this observation is that the summation over h in the 
partition function will be dominated by configurations near local minima due to a steepest 
descent phenomenon. For the rest of this section we shall assume that we are in a non- 
degenerate situation where the global minimum is unique and can be identified with the 
average height function (h) (we shall have to deal with degenerate situations later on). 

In order to go on we need an expression for the contribution to the free energy of 
configurations near a given height function h(x, y). We come to the main hypothesis of this 
reasoning: we assume that the only effect of the FBC is to create local polarizations. Thus, 
in every small patch V of the domain T> in which one can consider the local polarization 
P to be approximately constant, the free energy in T>' is a function of P only. Now we 
already encountered a model which is explicitly translationally invariant (i.e. with constant 
local polarization): the model with PBC, in which one can constrain the polarization to 
the same value P by Legendre transform. We therefore conjecture the equality of free 
energies per vertex in these two situations as we take the continuum limit. There are 
many consistency checks of this hypothesis, some of which we shall see below. The free 
energy being extensive, we reach the conclusion that the contribution of configurations 
near h(x,y) is given by 



defined earlier. Note again that since h is assumed to be 1-Lipschitz it is almost everywhere 
differentiable and the integral has a meaning. 




(2.9) 



where G(P X , P y ) = G{P) is the free energy per vertex with PBC at fixed polarization P 



Eq. ( |2.9| ) is all we need to provide our variational principle. It is clear now that 
the typical height function (h) is simply obtained by minimizing ( |2.9| ) over all 1-Lipschitz 
functions h satisfying the boundary condition h\gx> = ho. And the free energy per vertex 
is given by 




where A(T>) is the area of the domain T>. Eq. ( |2.10| ) is our key formula. 
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Minimizing the functional above leads to interesting questions. Indeed it is known 
that if G(P X , P y ) is a twice differentiable strictly convex function of P x and P y then h(x, y) 
is the unique solution of an elliptic partial differential equation (PDE), namely: 

d 2 G d 2 G d 2 G 

-Qp^{-d y h, d x h) d xx h - 2 — — (-d y h, d x h) d xy h + gp^i-dyh, d x h) d yy h = (2.11) 

It turns out that for no region of the parameters of the model is it fully the case; however 
the problems are more or less severe depending on the regime. The general principles of 
statistical mechanics ensure that the function G is convex, but there is the possibility of 
singularities (in our case, as we shall see, the endpoints \P X \ = \P y \ = 1 as well as the 
special point P x = P y = 0); and of linear parts (corresponding to the mixing of two phases 
at a first order transition) . We are thus led to the study of the phase diagram of the model, 
which will be the object of the next section. 

Let us conclude here by describing how the results for periodic and anti-periodic 
boundary conditions mentioned in the introduction can be recovered in a self-consistent 
way in our formalism. With PBC, one can still define a height function on a fundamental 
domain, say (x,y) e [0, L] x [0,1/]; however, the boundary conditions now read 

h(x = L, y) = h(x = 0, y) + LP y h(x,y = L') = h(x, y = 0) - L'P X (2.12) 

where the total polarization per vertex P is arbitrary. The free energy is now given by 
Eq. ( |2.10|) , but where h is subject to the boundary conditions (|2.12| ). Linear functions 
clearly satisfy these, and, for given P x and P y , minimize fl2.1U| ) due to the convexity of 
G. This means that the local polarization is constant and we immediately recover the 
free energy G(P); in particular the unconstrained free energy F, obtained by minimizing 
over P, is simply G(0). Modifiying slightly the argument leads to an identical result for 
anti-periodic or free boundary conditions. 



3. Bethe Ansatz solution and phase diagram 

The six-vertex model being integrable, one might think that a closed expression exists 
for the free energy F(E) of the six-vertex model with PBC in an electric field. In fact, 
even though Bethe Ansatz equations exist for an arbitrary field E, only for E = can they 
be solved in the thermodynamic limit. We briefly review here the Bethe Ansatz equations. 
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Consider the six vertex model with PBC on a rectangle of size Af = N x M. One 
always starts by using a transfer matrix formulation of the PBC partition function: 



Z = txT M (3.1) 

where T is the usual transfer matrix, acting on [<C 2 )® N that is on rows of vertical arrows 
of the lattice. The goal is to diagonalize T; the largest eigenvalue will then provide us with 
the free energy in the thermodynamic limit. Due to the conservation of arrows (and the 
PBC), the number n of up arrows is independent of the row and we can fix it. n being 
related to the vertical polarization by P y = 1 — 2n/N, the vertical eletric field plays no role 
and can be set to zero; so that the resulting free energy F = — limjv,M->-oo jjjr l°g^ is a 
mixed function of P y and E x , related by Legendre transform to both functions F(E X , E y ) 
and G(P x ,P y ) defined earlier. 

There are many equivalent diagonalization procedures. In the coordinate approach 
(see for example ||14j| ) , one makes an Ansatz on the wave function of the n up arrows. It 
is characterized by n distinct momenta kf, these have to satisfy constraints (Bethe Ansatz 
Equations, BAE) obtained by having one spin up circle around the strip of width N, 
interacting with all the other spins up: 



jk~N 



Y[B{kj,ki) j = l,...,n (3.2) 



e=i 

1*3 



where 



and 



1 + exp(4E x + ijp + g)) - 2A exp(2E x + ip) ( . 

[P:Q) l + exp(4E x + i(p + q))-2Aexp(2E x + iq) { ' } 



a 2 + b 2 - c 2 



2ab 

The eigenvalue of the transfer matrix is then given by 



(3.4) 



where 



T=(ae E *)"l[L(k j ) + (be- E *)"l[M(k j ) (3.5) 

3=1 3=1 



= a/6 + exp (a + 2 £ )-2A 

exp(ik + 2E x )a/b — 1 
= b /a + eM-ik-2E)-2A 
K ' exp(-ik-2E x )b/a-l K ' 
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In order to write down the final formula for the free energy, it is convenient to introduce 
some more notations. We set 



a = e £a b = e £b c = e £c 



(3.7) 



and split the contribution to the polarization of an arrow into its two endpoints, so that 
the energies of the six vertex configurations become 



ei = e a — E x — E y — >- 



A 



->- 



A 



£3 = e& — E x + E v — >- 



V 



->- 



V 



£5 = e c — 



->- 



A 



v 



■^T- 62 = £a + E x + E y 



A 



-<— €4 = €b + E X — Ey 



■<- ee = e c 



Then as M —> 00 one can write 



F = lim min 



1 n 1 n 

^^log|L(%)i,6 4 --^log|M(fc 



J=l 



(3.8) 



where n = -^-(l — P y ), and the kj are solutions of the Bethe Ansatz Eqs. ( |Q| ) chosen to 
minimize F. 

There is an obvious symmetry P y <-> —P y , and in what follows we shall always restrict 
ourselves to P y > 0, that is < n < N/2. This upper bound on n appears naturally in 
the BAE. Several scenarios can occur: the largest eigenvalue of T can be attained by the 
trivial state with all spins up (n = 0); it can correspond to the number n reaching its 
maximum value n = N/2; or it can correspond to a state with an intermediate number 
< n < N/2. After Legendre transformation, similar phenomena take place for P x . We 
are therefore led to the following definitions: 

o There may be regions of E for which the polarization is constant and satisfies \P X \ = 
\P y \ = 1. This corresponds to a frozen phase with all the arrows aligned: we call it the 
ferroelectric phase (F). We can distinguish which vertex configuration (from 1 to 4) is 
favored by the subscript: Fi, F 2 , F 3 , F4. 
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o There may be regions of E for which the polarization is constant and is zero: we call it 
the antiferroelectric phase (AF). 

o Finally, in regions of E where the polarization is non-constant, there is no particular 
order and we call this the disordered phase (D). 

The phase diagram depends very much on the value of the parameter A (defined by 
(|3.4j) ), and we shall examine three cases separately. Note that in what follows, we shall 
always assume a > b due to the symmetry between a and b. More detailed information on 
the phase diagram can be found in fl4[] . 

3.1. The regime A > 1 

This is one of the two "low temperature regimes". As in all other regimes, it seems 
impossible to derive exact expressions for the free energy in arbitrary field. Instead, it is 
natural to look for a low temperature expansion. Let us start at zero temperature, that is 
A = +00. In order to draw the zero temperature phase diagram, it is sufficient to compare 
the four energies and pick the lowest. This leads to Fig. 2. 



Fig. 2: Phase diagram as A — > +00. 

There is no bulk entropy and the free energy per vertex is simply F = min(ei, 62, £3, £4). 
The Legendre transformation is highly degenerate since the whole regions Fj , % = 1 , . . . , 4 
correspond to the points \P X \ = \P y \ = 1. The square (P x ,P y ) G [—1,1] X [— 1,1] is now 
divided into two regions, P x — P y ^ 0,@ which correspond to the two triple points of Fig. 2 
and where the phases Fi, F2, and F3 or F4 mix. In these regions G is linear, cf Fig. 3. 

2 Had we chosen b > a it would have been P x + P ' y ^ 0. 



E, 



y 




t b- t a 
2 
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F 2 

Fig. 3: Free energy as a function of the polarization at A = +00. 
The expression of the free energy: 

G(P X , P y ) =e a + 6j h^\P x - P y \ (3.9) 
leads, for FBC, to the minimization of the following quantity: 



/ = J J dxdy \d y h + d x h\ 

where we recall that h is defined on a convex domain V, is fixed to be ho on the boundary 
of £>, and satisfies \d x h\ < 1, \d y h\ < 1. 

We make the change of variables: u = x — y, v = x + y so that I = JJ dudv \d u h\. The 
integral over u is clearly minimum when d u h has fixed sign. Let us call ui(v) < u 2 {v) the 
intersections of dV and the v = est line; we see that for fixed boundary values h (ui(v)) 
and ho(u 2 (v)), the minimum is 

dv\h (u 2 (v)) - ho(u!(v))\ 



How many functions realize this minimum? One can always choose h to be linear in u at 
fixed v: 

hi \ u-ui(v) u-u 2 {v) , , 

h{u,v) = — —ho{u 2 {v)) — —hoimiv)) 

u 2 {v)-ui(v) ui(v)-u 2 (v) 

This function is 1-Lipschitz if h is. However, for a generic function ho, this will not be 
the only choice, and in fact there will be an infinity of possible h; this is related to the 
non-strict convexity of G, which is itself related to the fact that the phase transitions of 
Fig. 2 are first order. 

The final expression for the FBC free energy for FBC at A = +00 is: 

G FBC = e a + I -^-\ho{u 2 {v)) - ho{u x {v))\ (3.10) 



A(V) 

It is however important to understand that the degeneracy in the choice of h is an artefact 
of zero temperature, and is lifted by any finite value of A. Indeed, let us now look at a 
typical phase diagram when 1 < A < +00 (Fig. 4). 



10 



F 4 / 
——^^ \ 


y 


Fi 








- 


F 2 










1 Fa 



Fig. 4: Phase diagram in the A > 1 regime. 

There is a first order transition left between phases Fi and F 2 on a finite interval; 
however the triple points have been replaced with disordered regions with second order 
phase transitions at the boundary. The result, at the level of the Legendre transform, 
is that the two flat faces of G(P Xl P y ) acquire some curvature which make them strictly 
convex in the two regions P x — P y s$ 0. There remains only a one-dimensional flat valley 
between the two. Thus, we are led to the conclusion that in each of the spatial regions 
which correspond to P x — P y «s (which are determined, as stated earlier, by the sign of 
ho(u2(v)) — ho(ui(v))), h is either "saturated" (\d x h\ = 1 or \d y h\ = 1), or the solution of an 
elliptic PDE. The zero temperature limit of these equations is given by the first correction 
to the free energy. Here we shall not attempt the calculation of this correction since it 
is fairly involved, and a similar (though technically easier) calculation will be performed 
in the more interesting regime A < —1. Note finally, that in the (non-generic) case of 
boundary conditions such that there are regions for which ho(u2(v)) = ho(ui(v)), there is 
of course a unique way of minimizing the free energy for all A > 1, which corresponds to 
remaining in the flat valley (d x + d y )h = (mixture of Fi and F 2 ). 

3.2. The regime A < -1 

This is traditionally called the antiferroelectric regime, but in an electric field a dis- 
ordered phase and even ferroelectric phases can occur. 

We start once again with the zero temperature limit, which corresponds to A — > — oo. 
In this limit, and in zero field, the summation over six vertex configurations is dominated 
by the contribution of the ground state, which is, up to a reversal of arrows, given by 
Fig. 5, that is purely made of type c configurations. 
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Fig. 5: Antiferroelectric ground state. 
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Fig. 6: Phase diagram as A — > — oo. 

The full phase diagram at T = is simply obtained, just as in the case A — > +oo, by 
minimizing the energy of the six types of vertices, see Fig. 6. 

We note that there are first order phase transitions from antiferroelectric to ferroelec- 
tric phases. As before, we expect them to create degeneracies in the variational problem. 
In terms of the polarization P (see Fig. 7), we have: 

G(P X , P y ) =e c + e -^^\P x + P y \ + ejL ^\P, - Py\ (3-11) 



The domain V will now be divided into four subdomains depending on the sign of 
(d x — d y )h and (d x + d y )h. In each of these subdomains, the free energy being linear, there 
is generically an infinite number of functions h minimizing it and satisfying the constraints 
\d x h\ < 1, \d y h\ < 1. The explicit determination of the subdomains is however not as 
simple as in the ferroelectric phase and can only be achieved on a case-per-case basis. 
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Fig. 7: Free energy as a function of the polarization at A = — oo. 




Fig. 8: Phase diagram of the A < — 1 regime. 

Let us now switch on the temperature. When — oo < A < —1, as we can see on 
Fig. 8, there is an intermediate disordered phase which smooths the transitions to second 
order. 

It is convenient at this point to reintroduce a temperature T = 1/(3, so that the 
weights become 

a = e ~^/T h = e ~e b /T c = e - £c /T ( 3 _ 12 ) 

Similarly to the regime A > 1, we expect the first low temperature correction to lift the 
degeneracy and produce the PDE satisfied by h as T — > in the disordered regions. Let 
us now turn to this calculation. 

Because of the symmetries of the phase diagram 6, we can choose to expand around 
any of the four triple points. Note that these points are not equivalent from the point 
of view of the Bethe Ansatz. We choose E x near £a ~ £b ; as a consistency check, if P y is 
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allowed to vary one should recover a value of E y near <Ea ^ <Eb — e c . Also, this corresponds 
to the region \P X \ < P y (cf Fig. 7). 

First we consider the simplification of the Bethe Ansatz in the limit T — > 0. Assuming 
the exp(ikj) to remain bounded in this limit, we immediately obtain from Eq. ( |3.3| ) that 

B(p, q)~- e Kp-Q) 

so that the Bethe Ansatz Eqs. Q3.2|) have a straightforward solution. Explicitly, taking the 
logarithm and introducing half-integers Ij (so that the 21 j have the opposite parity of n), 
we find 

kj(N-n) = 2irIj-K 

where K = kj is the total momentum, which is itself quantized: K = 2tcI/N with 
I = ^2 a Ij , so that 

kj = _?^(J i -I/tf) (3.13) 



The ground state is obtained by choosing the kj as close to zero as possible, that is 

In particular the edge of the "Fermi sea" is ±Q where 

1 — P 

= 7T 



1 + Py 



Next we consider the behavior of the energy functions ( |3.6| ) as T — > 0. We introduce 
a scaling variable u such that 2i? x = e a — ej, + Tw; then Eqs. (|3.6| ) become 



c 2 /a6 
e u e lk — 1 



If it > ae Ex / T ^> be~ Ex / T so that it is the first term in Eq. ( |3.5|) , involving L(k), that 
dominates; whereas if u < it is the second term, involving M(k). 

Let us start with u > 0. Using Eqs. ( |3.8|) , ( |3. 14j) and (|3.15|) we find 

n f +Q dk , c 2 /a6 
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This is a dilogarithm integral; expanding in powers of e u results in 



oc 



F = ^P y (e a + e b -Tu) + (l- P y )e c -T^-^J2^^ sin ( m Q) u>0 ( 3 - 17 ) 



2tx ^— ' m 2 

Note that E y = = ^± - e c + 0(T), as expected. 

If u < 0, expanding M in powers of e u , one finds the same expression but with u 
replaced with —u. It is a nice consistency check of our approach that these two expressions 
are analytic continuations of each other as u crosses zero. Indeed it is clearly continuous 
at u = 0, and the two explicit expressions of its derivative 

2 dF { P y - i±5c argfl - e - u -^) u>0 



Pr = 



7T 



T du y -P y + i±parg(l - e u ~ lQ ) u<0 

turn out to be equal. In fact, the relation between u and P can be rewritten much more 
conveniently as 



SHI I 



sin I 7T \ + 



P X \<Py (3.18) 



When u — > +oo, P x ~ P y , and taking the Legendre transform of ( 3.16 ) we obtain 



G = P y e a + (1 - P y )e c P x =P y (3.19) 

which is what is expected (first order phase transition between phases Fi and AF). A 
similar results holds for u — > — oo. 

Finally, noting that u = Cb ^ Ca + 2-§p-, we see that Eq. (|3.18| ), supplemented by 



boundary condition ( |3.19| ), defines entirely G. The explicit expression involves once again 
dilogarithms and will not be written down here. What is more interesting is the second 
derivatives of G, which provide the PDE ( |2.11| ) satisfied by h. After a few calculations we 
find: 

2ttP v 



d 2 G n sin 



8P* 2 (1 + Py)sin ^) sin ^) 
d 2 G n sin |^ - P x sin * 



dP x dP y 2 (i + Py) 2 sin H^+EA sin 



2 " 2 tv(1-P x 



( ■ 2-irPr t-> ■ 2tt \ , • 2 ir(l+P x ) ■ 

d 2 G _ n [ sm j+jg ~ P * sm 1+7^ J + sm ^+7^ sm i+P, 

dP 2 ~ 2 M I P V^.in 27TP y "in l(l±^) „: n ^(1--PQ 

w U + ^yJ sm i+p h sm 1+Py sm i+p h 



(3.20a) 
(3.206) 

(3.20c) 
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One has f^fpl - y dP x dP y ) = (i+p v ) 4 ; ^ e corres P on ding PDE, of the form ( |2.11|) , 
is elliptic. 

Note that Eq. ( |3.20| ) are only valid for \P X \ < P y ; similar expressions can be found 
in the other three regions by using the various discrete symmetries of the model. These 
expressions become singular when one approaches the boundaries \P X \ = \P y \; this is to be 
expected in a perturbative treatment around T = 0, but such singularities should disappear 
non-perturbatively. On the contrary, the singularity at P x = P y = should remain due to 
the nature of the phase diagram (existence of an antiferroelectric phase, cf Fig. 8). 

3.3. The regime -1 < A < 1 

This is traditionally called the disordered regime; however, as can be seen on the phase 
diagram (Fig. 9) , once the electric field is turned on both disordered and ferroelectric phases 
occur. 




Fig. 9: Phase diagram for — 1 < A < 1. 

The only difference with the phase diagram for A > 1 is that there is no more phase 
boundary at E x = E y = 0; instead the neighborhood of this point is entirely inside the 
disordered phase. The transitions from disordered to ferroelectric phases are second order. 
This means that we expect that for given FBC, the domain T> will be divided into several 
zones (with definite boundaries), each corresponding to a phase of Fig. 9: typically the 
ferroelectric phases will lie close to the boundaries with maximum slope (|<9 x /io| = 1 or 
\d y ho |=1), where their effect is strongest, while the disordered phase occupies the remaining 
space. Inside the disordered zone, the function h will satisfy an elliptic PDE. 
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Unfortunately, as in other phases this cannot be made explicit except in one particular 
case: the "free fermion point" , which will be described separately in the next section. For 
a general value of A between —1 and 1 the only exact result that seems derivable is the low 
field expansion of the free energy. This is a standard computation which will be sketched 
now. 

First we introduce another parameterization of the weights: 

a = sin^-(l — C) 6 = sin^-(l + C) c = sin7 (3-21) 

with — 1 < C < 1 (or — 1 < C < once we restrict ourselves to a > 6), < 7 < tt, so that 
A = — cos 7, and make the standard change of variables in Eqs. (|3.2[) -( 



Ak 3 +2E m _ sinh 7 (Aj+z/2) 

6 " sinh 7 (A J -z/2) [6 - ZZ) 



The BAE become: 



smh-f(X j + i/2) c _ 2Ex 
sinh 7 (Aj — i/2) 



N - TT smh 7 (A J -A, + z) 

-llsinh7(A,-A,-z) ^- 1 ""-' n ^-23) 

1^3 



and the functions L and M in terms of which the free energy is express are given by 

Note that in the rhs of Eqs. ( |3.23j ), only the differences Xj — Xi appear. 

Let us now assume E x = 0. Then equations ( |3.2^j| ) are real and their "ground state" 
solution is given by real A^. If we further assume that n = N/2 (that is after Legendre 
transformation, E y = 0), then it is known that the A^ fill the entire real axis in the limit 
iV — > 00, and their density is given by a simple integral equation. Namely, by taking 
the logarithm of ( |3.2ij| ) and differentiating, we find that the density po(X), normalized by 
J po(X)d\ = n/N, satisfies 

K* Po = (j) (3.25) 

where ★ denotes convolution product, and we give K and <p directly in Fourier transform: 
K{k) = fK(X)eM^X)dX = 2 sinh ^- i ;; ) r h - and 0(«) = so that we 

simply find p («) = or 

Po(A) = - 1 (3.26) 
2 cosh7rA 
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If n < N/2, the equation remains similar, but only an interval of the real axis, of the 
form [—Q,Q], is filled. The integral equation becomes equivalent to a 2 x 2 Riemann- 
Hilbert problem which is in general unsolvable. However in the limit Q — > oo (that is P y , 
E y — > 0), the effect of the interval endpoints on each other become negligible are we are left 
with a simple Wiener-Hopf problem. Let us briefly recall the principle. We concentrate 
on one of the two edges, say A near Q, consider the new density p(X) and notice that its 
Fourier transform p(n) exp(—2inQ) is holomorphic in the half-plane Imn < and bounded 
by a polynomial at infinity. The equation in Fourier transform is now: 

K-\ K )a(K) + p(n) =pq{k) (3.27) 

where the unknown function cr(A) (the density of "holes") is, in our approximation, such 
that <x(«) exp(— 2ikQ) is on the contrary holomorphic in the half-plane Im/t > 0. Solving 
this equation involves the decomposition K(k) = K+(k) / K-(k) where 

K+(k) = — — T( : K {l] , r^e aK (3.28a) 

+ V ' r((l/7- l/7r)m)r(l/2 + m/yr) v ; 

,_ (k) = , r (1 -(i/,-i/^ : )r( 1 /2-^) e „. {3m) 

and a = — l/?r) log(l — ■y/Tv) log(7/-7r). In general one must also decompose the 
right hand side of ( |3.27|) , but since Q —>■ oo one can consider only the contribution of the 
dominant pole at k = —in/2, so that the solution of this analytic problem is elementary: 

fi( K ) - e 2iQ(K+iir/2) ^+(~W 2 ) i (o no) 

p(k) is undefined since in this Q — ► oo limit, p(X + Q) diverges exponentially as A — ► — oo; 
however one can write 

p( K ) = e 2 'QW) K + {-iv/2) -i 

K + {k) ^ + z7r/2| deformed 

where it is understood that to recover p(X) one must integrate on a deformed contour that 
goes around the pole k = —iir/2. 

To recover the expression of F as a function of P y we write, remembering that there 
is also a contribution around A = — Q: 



P y = 2 /(po(A)-p(A))dA=-2p(0) 



A Q K + {-i«/2) 
K+(0) 
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and 

F = F + 2Re J (p (A) - p(X)) logL(A)dA 

= F + 2-Res« =i7r/2/ 5(K) fo£~L(-m/2) + ■■■ (3.32) 

7T 

so that after some calculations we find our first correction to the free energy: 

Fi = J(1-7A)oos(Ctt/2)^ (3.33) 

If we now switch on the horizontal field E x , the factor exp(— 2E X ) in Eqs. ( |3.23| ) breaks 
their reality and we can no longer assume that the A^ are real; instead they will fill some 
curve in the complex plane. At first order in E x the Xi will be displaced by a purely 
imaginary value, so that we can write 

Xi = X° l + iE x u(\i) 

where the A^ are the solutions of the BAE at E x = 0. Expanding the BAE we easily find 
the equation for the function p u (X) = p(X)u(X); in Fourier transform, 

K(k)p u (k) + a u (n) = 28 (k) 

with the usual analyticity constraints on p u and a u . The solution is simply: 

a u (n) = ^ (3-346) 

This time the correction to the free energy is 

F 2 = -E x + 2iE x J &X Pu (X)(\ogL)\X) + 0(E 2 X ) (3.35) 

Rewriting the integral in Fourier transform, we see that its first two dominant contributions 
come from the poles at k = and k = ztt/2. The pole at k — compensates exactly the 
trivial term —E x , and we are left with: 



F 2 = E x P y sin(C7r/2) (3.36) 
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Note that the absolute normalization of this term is meaningful. 

The computation can be carried out to the next order in P y \ however it is simpler 
to obtain it by P x «-> P y symmetry. Combining Eqs. ( |3.33| ) and ( |3.36| ) and taking the 
Legendre transform, we find: 

G(P Xl P y ) = G(0, 0) + ^(1 - 7/tt) ^^2) [ P - + 2 sin (C^/2)^^ + Py] + ■ ■ ■ (3-37) 

This is the low polarization expansion of G. 

The corresponding linearized equation for h is: 

Q2 Q2 Q 2 \ 

2 sin(C7r/2)— — + )h = (3.38) 



<9x 2 dxdy dy' 

This is a very simple result, which has the following interpretation: for a low field we 
expect the usual conformal invar iance to hold in the continuum limit. It is known that the 
phase |A| < 1 is critical and its infrared limit is described by a free massless boson. Here, 
h plays the role of the euclidean bosonic field. Therefore it must satisfy Ah = 0. However 
it is well-known that the spectral parameter £ creates an anisotropy of the lattice: the two 
directions of the lattice must be considered as forming an angle of ^(1 — £), which leads 
to Eq. ( g3g ). 

3. 4- The free fermion point A = 

When A = 0, that is when the Boltmann weights satisfy the relation a 2 + b 2 = c 2 , the 
system describes free fermions; this is particularly clear when one considers the expression 
( p.3|) of the two-body interaction in the Bethe Ansatz equations, which becomes simply 
B = — 1; the BAE ( |3.2| ) themselves become 

e ikjN = (_!)»-! (3.39) 

i.e. describe (ignoring the sign issue related to our "bosonic" description of fermions) the 
quantization of the momenta of free particles. The kj still satisfy an exclusion principle, 
and the largest eigenvalue of the transfer matrix is obtained by choosing 

The edge of the Fermi sea is 

Q= 7 L{l- Py ) (3.41) 
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and the free energy is given by an integral of a type encountered above: 



Q 



Q 



dk 
2rv~ 



log 



a / b + e ik+2E x 

_e ik+2E *a/b-l 



(3.42) 



Finally, by taking the Legendre transform and differentiating twice, we obtain the 
PDE (of the type ( |2.11| )) satisfied by h in the unsaturated regions at A = 0: 







f 2 , i2\ 2 n ^x , 

[a + b ) cos — — h n 



2 n(h x + h y ) 2 ir(h x - h 
+ 2 [ a cos — - b cos 



yj 



TV 



h,, 



TV 



h. 



cos 



cos 



h 



xy 



+ (a 2 + b 2 ) cos' 



ivh 



h 

2 w 



(3.43) 



It is easy to check that this equation is elliptic for all \h x \ < 1, \h y \ < 1. If one sets a = b, 
after rotation of tv/4 one recovers the PDE found in for domino tilings. 



4. Conclusion 

We have proposed in this paper a variational principle (Eq. (|2.10|) ) for the the six- 
vertex model in the continuum limit. For generic fixed boundary conditions this variational 
principle has a non-trivial solution, so that the existence of a thermodynamic limit for 
the six-vertex model is invalidated. A qualitative discussion of the phase diagram in 
the various regimes of the model has suggested that in general the system undergoes 
phase separation into regions which belong to one of three phases: two ordered phases 
(ferroelectric, anti-ferroelectric), and one disordered phase in which one can write a partial 
differential equation satisfied by the height function (Eq. ( |2.11| )). This leads to a first 
open question: it would be nice to have a rigorous proof of the existence of such a phase 
separation. In particular, it is not obvious at all in our framework why an anti-ferroelectric 
region should exist at all, though we conjecture it is generically the case for A < — 1. 
The calculations presented in this paper are essentially a perturbative analysis around 
zero temperature, which is highly degenerate and does not help to solve this issue. We 
have however given the explicit PDE in the disordered regions in the limit A — » — oo 
(Eqs. (3.20)). 

It is unclear at the moment how useful this variational principle is in given particular 
cases. One application is the domain wall boundary conditions mentioned in the intro- 
duction. Note that for A = 0, a = b, this is equivalent to the problem of domino tilings 
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of the aztec diamond H, for which the height function is known exactly [11] and solves 
Eq. ( p.43|) . For A > 1, a > b, using the formalism of section 3.1, it is very easy to proof 
the conjecture of section 4.1 (Fig. 4) of 0. In our language, it states that the height 
function is given for all A > 1 by h(x,y) = \x + y\, — 1 < x,y < 1 (i.e. it minimizes the 
free energy among the functions with identical boundary values at \x\ = 1, \y\ = 1); it 
corresponds to the non-generic case mentioned at the end of section 3.1 where one remains 
in the flat one-dimensional valley P x = P y . Similarly, for A —>■ oo, one can prove the 
conjecture of section 5.3 (Fig. 3) of ||, with mild assumptions on the possible form of the 
AF/F domains. One should however note that the DWBC, with their maximum slope on 
all boundaries, are very non-generic boundary conditions, which explains in particular the 
absence of degeneracy for A = ±oo. Other similar boundary conditions which give rise 



to exact determinant formulae are introduced in [ 15(1 , and can probably be analyzed along 
the lines of @ or || ; however, they presumably give rise to the same height function, since 
they can be obtained from DWBC configurations by dividing by a certain discrete symme- 
try which is not spontaneously broken. It would therefore be useful to find other exactly 
solvable boundary conditions to provide non-trivial examples of the variational principle 
developed here. 
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